import unittest as UT
import numpy as np
from numpy import abs
from numpy.random import random
from Bernstein import *

j = complex(0,1)

class test_Bernstein(UT.TestCase):

    def setUp(self):
        self.rP = (random((1,2)).squeeze(), random((1,2)).squeeze(),
                   random((1,2)).squeeze(), random((1,2)).squeeze())
        self.tRange = np.linspace(0, 1, 11)
        
    def testLinear(self):
        checkLine = (1-self.tRange)*self.rP[0].reshape(-1,1) + \
                                    self.tRange*self.rP[1].reshape(-1,1)
        LB = (Bernstein([self.rP[1], self.rP[2]]),
              Bernstein([self.rP[2], self.rP[0]]),
              Bernstein([self.rP[0], self.rP[1]]))
       # Test s(t, i)
        self.assert_(LB[0].s(self.rP[1][0], self.rP[2][0]) == \
                          self.rP[1][0] - self.rP[2][0])
       # Test b, B
        self.assert_(LB[2].b(0.4, 0)==0.6*self.rP[0][0] + 0.4*self.rP[1][0])
        self.assert_(LB[2].b(0.6, 1)==0.4*self.rP[0][1] + 0.6*self.rP[1][1])
        self.assert_((LB[2].b(0).T==self.rP[0]).all())
        self.assert_((LB[2].b(1).T==self.rP[1]).all())
        self.assert_(np.isclose(LB[2].b(1.0/3).T,
                                        2*self.rP[0]/3 + self.rP[1]/3).all())
        self.assert_(np.isclose(LB[2].B(self.tRange), checkLine).all())
      # Test bprime & Bprime
        self.assert_(LB[2].bprime(random(),0) == \
                                             (self.rP[1][0] - self.rP[0][0]))
        self.assert_(LB[2].bprime(random(),1) == \
                                             (self.rP[1][1] - self.rP[0][1]))
        self.assert_((LB[2].bprime(random()).T == \
                                            (self.rP[1] - self.rP[0])).all())
        self.assert_((LB[2].Bprime(checkLine).T == \
                                            (self.rP[1] - self.rP[0])).all())
      # Test SB
##        a = 0
##        b = 1
##        z0 = 0 + 0*j
##        z1 = 1 + 0*j
##        A = (b**2 - a**2) * z1 / 2
##        B = (1 - (b + a)/2) * z0
##        tmp = LB[0].SB(0, 1, 0)
##        print tmp, (z1 + z0) / 2
##        self.assert_(np.isclose(tmp, A + B).all())
    
      # Test length
        for n in range(3):
            dz = np.dot((self.rP[(n+2)%3] - self.rP[(n+1)%3]).T, [1,j])
            self.assert_(np.isclose(LB[n].length(), abs(dz)))
      # Test tGivenCoord
        t = random()
        lb = LB[2].b(t)
        self.assert_(np.isclose(LB[2].tGivenCoord(lb[0], 0), t))
        self.assert_(np.isclose(LB[2].tGivenCoord(lb[1], 1), t))
        
      # Test whereConnected, isClosed
        openPath = Path([LB[2], LB[0], Bernstein([self.rP[2], self.rP[3]])])
        self.assert_((openPath.whereConnected() == 
                                        np.array((0,1,1), dtype=bool)).all())
        self.assert_(not openPath.isClosed())
        osedPath = Path([LB[0], LB[1], LB[2]])
        self.assert_((osedPath.whereConnected() ==
                                         np.array((1,1,1),dtype=bool)).all())
        self.assert_(osedPath.isClosed())
      # Test reparam, integrate, isInside, not isInside
##        start, end = osedPath()
##        self.assert_

        self.assert_(not openPath.isInside(random((1,2)).squeeze()))
        tol = 1e-12
        for n in range(3):
            testLine = Bernstein([eval('self.rP[' + str(n) + ']'),
                (1-t) * eval('self.rP[' + str((n+1)%3) + ']') + \
                t * eval('self.rP[' + str((n+2)%3) + ']')])
            testInP1 = testLine.b(tol).squeeze()
            testInP2 = testLine.b(1 - tol).squeeze()
##            self.assert_(osedPath.isInside(testInP1))
##            self.assert_(osedPath.isInside(testInP2))
            testOuP1 = testLine.b(-tol).squeeze()
            testOuP2 = testLine.b(1 + tol).squeeze()
##            self.assert_(not osedPath.isInside(testOuP1))
##            self.assert_(not osedPath.isInside(testOuP2))
      # Test reparam
        


    def testQuad(self):
        qP = ((1,0), (0,1), (-1,0), (0,-1))
        checkQuad = np.power(1-self.tRange,2) * self.rP[0].reshape(-1,1) + \
                    2*self.tRange*(1-self.tRange) * self.rP[1].reshape(-1,1)\
                    + np.power(self.tRange,2) * self.rP[2].reshape(-1,1)
        QB = (Bernstein([self.rP[0], self.rP[1], self.rP[2]]), \
              Bernstein([self.rP[1], self.rP[2], self.rP[0]]), \
              Bernstein([self.rP[2], self.rP[0], self.rP[1]]))
       # Test s(t, i)
        self.assert_(QB[0].s(qP[1][0], qP[3][0]) == \
                             qP[1][0] - qP[3][0])
       # Test b, B
        t = 0.36*self.rP[0][0] + 0.48*self.rP[1][0] + 0.16*self.rP[2][0]
        self.assert_(np.isclose(QB[0].b(0.4,0), t))
        t = 0.16*self.rP[0][1] + 0.48*self.rP[1][1] + 0.36*self.rP[2][1]
        self.assert_(np.isclose(QB[0].b(0.6,1), t))
        self.assert_((QB[0].b(0).T==self.rP[0]).all())
        self.assert_((QB[0].b(1).T==self.rP[2]).all())
        self.assert_(np.isclose(QB[0].b(1.0/3).T, 
            4*(self.rP[0]+self.rP[1])/9 + self.rP[2]/9).all())
        self.assert_(np.isclose(QB[0].B(self.tRange), checkQuad).all())
      # Test bprime & Bprime
        t = random()
        self.assert_(np.isclose(QB[0].bprime(t,0), 
                                       2 * (self.rP[1][0] - self.rP[0][0] + \
                    t * (self.rP[0][0] - 2*self.rP[1][0] + self.rP[2][0]))))
        self.assert_(np.isclose(QB[0].bprime(t,1),
                                       2 * (self.rP[1][1] - self.rP[0][1] + \
                    t * (self.rP[0][1] - 2*self.rP[1][1] + self.rP[2][1]))))
        self.assert_(np.isclose(QB[0].bprime(t).T, 
                                             2 * (self.rP[1] - self.rP[0] + \
                       t * (self.rP[0] - 2*self.rP[1] + self.rP[2]))).all())
        checkDeriv = 2 * ((self.rP[1] - self.rP[0]).reshape(-1,1) + \
          self.tRange * (self.rP[0]-2*self.rP[1]+self.rP[2]).reshape(-1,1))
        self.assert_(np.isclose(QB[0].Bprime(self.tRange), checkDeriv).all())
      # Test length
        rt2 = np.sqrt(2)
        testLength = rt2 + np.log(rt2 + 1)
        for n in range(4):
            tmp = Bernstein([qP[n], qP[(n+1)%4], qP[(n+2)%4]])
            self.assert_(np.isclose(tmp.length(), testLength))
      # Test tGivenCoord
        for n in range(len(self.rP[0])):
            x0 = self.rP[0][n]
            x1 = self.rP[1][n]
            x2 = self.rP[2][n]
           # First test for -discriminant case
            numer = x0*x2 - x1*x1
            denom = x0 - 2*x1 + x2
            if denom > 0: # test x < numer/denom to assure a neg. disc.
                x = numer/denom - 1
            elif denom < 0: # test x > numer/denom to assure a neg. disc.
                x = numer/denom + 1
            elif numer < 0: # if denom = 0 and numer < 0
                print "numer < 0: ", numer
                x = random() # any value of x should guarantee a neg. disc.
            else: # if denom = 0 and numer >= 0
                print "denom = 0 & numer >= 0: ", denom, numer
                x = np.NAN # no x will give a neg. disc => skip test
            if np.isfinite(x):
                t = QB[0].tGivenCoord(x, n)
                self.assert_(np.isnan(t))
       # Now test non-degenerate case
        t = random()
        qb = QB[0].b(t)
        test0 = QB[0].tGivenCoord(qb[0], 0)
        if len(test0)==1 and np.isfinite(test0):
            self.assert_(np.isclose(test0, t))
        elif len(test0) > 1:
            self.assert_(np.isclose(test0, t).any())
        else:
            self.assert_(np.isnan(test0))
        test1 = QB[0].tGivenCoord(qb[1], 1)
        if len(test1)==1 and np.isfinite(test1):
            self.assert_(np.isclose(test1, t))
        elif len(test1) > 1:
            self.assert_(np.isclose(test1, t).any())
        else:
            self.assert_(np.isnan(test1))
      # Test whereConnected, isClosed, isInside, not isInside
        openPath = Path([QB[2], QB[1], 
            Bernstein([QB[0]._P[0], self.rP[1], self.rP[3]])])
        self.assert_((openPath.whereConnected() == 
                                        np.array((0,1,1), dtype=bool)).all())
        self.assert_(not openPath.isClosed())
        osedPath = Path([QB[2], QB[1], QB[0]])
        self.assert_((osedPath.whereConnected() ==
                                         np.array((1,1,1),dtype=bool)).all())
        self.assert_(osedPath.isClosed())

    def testCube(self):
        cP = ((0,0), (1,0), (1,1), (0,1))
       # Test b, B
        checkCube = np.power(1-self.tRange,3) * self.rP[0].reshape(-1,1) + \
                              3 * self.tRange * np.power(1-self.tRange,2) * \
                  self.rP[1].reshape(-1,1) + 3 * np.power(self.tRange,2) * \
                              (1-self.tRange) * self.rP[2].reshape(-1,1) + \
                           np.power(self.tRange,3) * self.rP[3].reshape(-1,1)
        CB = (Bernstein([self.rP[0], self.rP[1], self.rP[2], self.rP[3]]),
              Bernstein([self.rP[3], self.rP[0], self.rP[1], self.rP[2]]),
              Bernstein([self.rP[2], self.rP[3], self.rP[1], self.rP[0]]))
       # Test s(t, i)
        self.assert_(CB[0].s(cP[1][0], cP[3][0]) == \
                          cP[1][0] - cP[3][0])
        x004 = 0.216 * self.rP[0][0] + 0.432 * self.rP[1][0] + \
               0.288 * self.rP[2][0] + 0.064 * self.rP[3][0]
        self.assert_(np.isclose(CB[0].b(0.4, 0), x004))
        x106 = 0.064 * self.rP[0][1] + 0.288 * self.rP[1][1] + \
               0.432 * self.rP[2][1] + 0.216 * self.rP[3][1]
        self.assert_(np.isclose(CB[0].b(0.6, 1), x106))
        self.assert_((CB[0].b(0).T==self.rP[0]).all())
        self.assert_((CB[0].b(1).T==self.rP[3]).all())
        self.assert_(np.isclose(CB[0].b(0.6).T, 
            (8*self.rP[0] + 36*self.rP[1] + 54*self.rP[2] + 27*self.rP[3]) \
                                                                 / 125).all())
        self.assert_(np.isclose(CB[0].B(self.tRange), checkCube).all())
      # Test bprime & Bprime 
        t = random()
        self.assert_(np.isclose(CB[0].bprime(t,0),
          3 * (self.rP[1][0] - self.rP[0][0] + \
          2 * t * (self.rP[0][0] - 2 * self.rP[1][0] + self.rP[2][0]) + t*t\
          *(self.rP[3][0]-self.rP[0][0]+3*(self.rP[1][0]-self.rP[2][0])))))
        self.assert_(np.isclose(CB[0].bprime(t,1),
          3 * (self.rP[1][1] - self.rP[0][1] + \
          2 * t * (self.rP[0][1] - 2 * self.rP[1][1] + self.rP[2][1]) + t*t\
          *(self.rP[3][1]-self.rP[0][1]+3*(self.rP[1][1]-self.rP[2][1])))))
        self.assert_(np.isclose(CB[0].bprime(t).T,
          3 * (self.rP[1] - self.rP[0] + \
          2 * t * (self.rP[0] - 2 * self.rP[1] + self.rP[2]) + t * t * \
          (self.rP[3] - self.rP[0] + 3*(self.rP[1] - self.rP[2])))).all())
        checkDeriv = 3 * ((self.rP[1] - self.rP[0]).reshape(-1,1) + 2 * \
          self.tRange * (self.rP[0]-2*self.rP[1]+self.rP[2]).reshape(-1,1) \
                   + self.tRange * self.tRange * (self.rP[3] - self.rP[0] + \
                                3 * (self.rP[1] - self.rP[2])).reshape(-1,1))
        self.assert_(np.isclose(CB[0].Bprime(self.tRange), checkDeriv).all())
      # Test length
        for n in range(4):
            tmp = Bernstein([cP[n], cP[(n+1)%4], cP[(n+2)%4], cP[(n+3)%4]])
            self.assert_(np.isclose(tmp.length(), 2))
      # Test tGivenCoord
        cb = CB[0].b(t).squeeze()
        print "t_in: ", t
        print "cb: ", cb
        test = CB[0].tGivenCoord(cb[0], 0)
        self.assert_((len(test)==1) and (np.abs(test-0.5) <= 0.5))
        test = CB[0].tGivenCoord(cb[1], 1)
        self.assert_((len(test)==1) and (np.abs(test-0.5) <= 0.5))
        

      # Test whereConnected, isClosed, isInside, not isInside
        openPath = Path([CB[0], CB[1], 
            Bernstein([self.rP[2], self.rP[3], self.rP[0], self.rP[1]])])
        self.assert_((openPath.whereConnected() == 
                                        np.array((0,1,1), dtype=bool)).all())
        self.assert_(not openPath.isClosed())
        osedPath = Path([CB[0], CB[1], CB[2]])
        self.assert_((osedPath.whereConnected() ==
                                         np.array((1,1,1),dtype=bool)).all())
        self.assert_(osedPath.isClosed())

    def tearDown(self):
        pass
        
##class test_Path(UT.TestCase):
##    
##    def setUp(self):
##        self.P0 = np.array([1,0])
##        self.P1 = np.array([0,1])
###        self.P2 = -self.P0
##        self.pathOpen = Path([Bernstein([self.P0, self.P1]), 
##                              Bernstein([self.P1, -self.P0])])
##        self.pathClosed = Path([Bernstein([self.P0, self.P1, -self.P0]), 
##                                Bernstein([-self.P0, -self.P1, self.P0])])
##
##    def test_isInside(self):
##        interiorP = [(0,0), (1 - 1e-12, 0), (1e-12 - 1, 0),
##                     (0, 0.5-1e-12), (0, 1e-12-0.5),
##                     (0.5, 0.125-1e-12), (-0.5, 1e-12-0.125)]
##        for P in interiorP:
##            self.assert_(self.pathClosed.isInside(P))
##
##        exteriorP = [(1 + 1e-12, 0), (-1 - 1e-12, 0),
##                     (0, 0.5 + 1e-12), (0, -0.5 - 1e-12),
##                     (0.5, 0.125 + 1e-12), (-0.5, -0.125 - 1e-12)]
##        for P in exteriorP:
##            self.assert_(not self.pathClosed.isInside(P))
        

##    def test_isOn(self):
##        self.assert_(self.pathOpen.isOn(self.P0+1e-7))
##        self.assert_(not self.pathOpen.isOn(self.P0+1e-5))
##        self.assert_(self.pathOpen.isOn(np.array((0.5,0.5))))
        

##    def tearDown(self):
##        pass
    
if __name__ == '__main__':
    UT.main()
